How 'sticky' are short-range square-well fluids? 
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The aim of this work is to investigate to what extent the structural properties of a short-range 
square-well (SW) fluid of range A at a given packing fraction n and reduced temperature T* = ksT/e 
can be represented by those of a sticky-hard-sphere (SHS) fluid at the same packing fraction and an 
effective stickiness parameter t(T*,A). Such an equivalence cannot hold for the radial distribution 
function g(r) since this function has a delta singularity at contact (r = ct) in the SHS case, while 
it has a jump discontinuity at r = Act in the SW case. Therefore, the equivalence is explored with 
the cavity function y(r), i.e., we assume that ysw{r\rj, T*; A) ~ ysns{r\r], t(T* , A)). Optimization 
of the agreement between t/sw and j/shs to first order in density suggests the choice r(T*,A) = 
[12(e 1/T * - 1)(A - We have performed Monte Carlo (MC) simulations of the SW fluid for 

A = 1.05, 1.02, and 1.01 at several densities and temperatures T* such that r(T* , A) = 0.13, 0.2, and 
0.5. The resulting cavity functions have been compared with MC data of SHS fluids obtained by 
Miller and Frenkel [J. Phys: Cond. Matter 16, S4901 (2004)]. Although, at given values of n and r, 
some local discrepancies between ysw and ysHS exist (especially for A = 1.05), the SW data converge 
smoothly toward the SHS values as A — 1 decreases. In fact, precursors of the singularities of ysHS 
at certain distances due to geometrical arrangements are clearly observed in ysw- The approximate 
mapping ysw — ► ?/shs is exploited to estimate the internal energy and structure factor of the SW 
fluid from those of the SHS fluid. Taking for j/shs the solution of the Percus-Yevick equation as 
well as the rational- function approximation, the radial distribution function g(r) of the SW fluid is 
theoretically estimated and a good agreement with our MC simulations is found. Finally, a similar 
study is carried out for short-range SW fluid mixtures. 



I. INTRODUCTION 

It is well known that colloidal particles in a sus- 
pension of free polymers interact through effective at- 
tractive forces (of entropic origin) with a range and 
strength determined by the size and concentration of the 
polymersi^ A simple model describing this effective in- 
teraction is the square-well (SW) potential 



co, r < a, 
Vsw(0 = ^ -e, a <r < Act, 
0, r > Act, 



(1.1) 



which accounts for excluded volume effects associated 
with the hard-core diameter ct plus an attractive layer 
of relative width A — 1 and strength e. Although the 
SW potential (with A ~ 1.5) was originally introduced as 
a simple model for normal liquids^ it has become even 
more useful in colloidal systems, where one typically has 
A<1.1. 

The interaction parameters ct and e can be used to 
define the length and energy units, respectively. Thus, 
the concentration of particles can be characterized by 
the volume fraction 77 = -|/Oct 3 , where p is the number 
density, and the temperature T can be measured in units 
of e as T* = ksT/e, where Ub is the Boltzmann constant. 
Therefore, at given r\ and T*, only the range A remains 
as a free parameter. This degree of freedom disappears 
in the so-called sticky-hard-sphere (SHS) limit^ where 



one takes the combined limits e — > 00 (i.e., T* — > 0) and 
A — > 1, while keeping constant the stickiness parameter 



12e 



l/T" 



(A-l). 



(1.2) 



The parameter r plays in the SHS fluid a role equiva- 
lent to that played by the reduced temperature T* in 
the SW fluid. Although, strictly speaking, a monodis- 
perse system of SHS is not thermodynamically sta- 
ble, a small degree of polydispersity is sufficient to re- 
store stability. 5 Given its simplicity and the fact that 
it can be exactly solved within the Percus-Yevick (PY) 
approximation^^ the SHS model has received con- 
siderable attention as a convenient model of colloidal 
suspensions 7 - 8 - 9 - 10 - 11 - 12 - 13 - 14 - 15 - 16 - 17 - 18 - 19 - 20 - 21 - 22 - 23 - 24 

It is obvious that the thermodynamic and structural 
properties of the SW fluid must approach those of the 
SHS fluid as A — > 1 at fixed 77 and r. This has been 
tested, for instance, in the case of the vapour-liquid crit- 
ical point An interesting problem not yet addressed in 
detail is the rate of the convergence from SW to SHS. In 
other words, how small must the width A — 1 be for the 
properties of the SW system to be practically the same as 
those of an equivalent SHS system? Comparison between 
Monte Carlo (MC) simulations^ for an SW fluid with 
A = 1.1 and theoretical results for SHS fluids shows a 
good agreement between the respective structure factors 
S(k)r2£££l although a systematic phase shift between 
both structure factors exists. This is a reflection of the 
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fact that the radial distribution function g(r), which is 
essentially the inverse Fourier transform of S(k), has im- 
portant qualitative differences in both classes of systems. 
While in the case of an SW fluid, g{r) is finite at r = ct + 
and presents a jump discontinuity at r = Act, a delta- 
peak singularity at r — ct + , followed by other delta and 
jump singularities at certain characteristic distanceSf&Si 
appear in the SHS case. 

The aim of this paper is two-fold. First, we want to 
assess to what degree a short-range SW fluid behaves 
as a suitably chosen "equivalent" SHS fluid. To that 
end, we have performed Monte Carlo (MC) simulations 
of SW systems with A = f .05, A = 1.02, and A = 1.01, 
and have compared their structural properties with re- 
cent MC simulations 21 for SHS fluids at the same pack- 
ing fraction and effective stickiness. The comparison is 
not carried out with the radial distribution function g(r), 
but with the much more regular cavity (or background) 



function y(r) 



oV (r)/k B T 



g(r). It is observed that, as 



expected, the SW cavity functions converge toward the 
SHS ones as A decreases. However, the results show that 
the range A = 1.05 cannot be considered small enough 
to get a cavity function hardly distinguishable from the 
SHS one. The second point we address is the estimate 
of the structure factor and thermodynamic properties of 
short-range SW fluids from the knowledge of the prop- 
erties of SHS fluids. As said before, the starting point 
is the assumption that both cavity functions are approx- 
imately the same. The approximate estimates for the 
internal energy and structure factor of the SW systems 
are tested against actual MC simulations, finding a good 
agreement. Also in this context, we take advantage of the 
exact solution of the PY equation for SHS^ as well as of 
a more refined theory^ to propose analytic approxima- 
tions for the radial distribution function of SW fluids. A 
parallel study is carried out for short-range SW fluid mix- 
tures. Due to the scarcity of simulation results for SHS 
mixtures, we compare our MC results with the solution 
of the PY theory. 6 It is observed that in the case of mix- 
tures the influence of A ^ 1 on the mapping SW— >SHS 
is similar to that in the one-component case. 

This paper is organized as follows. The basic equa- 
tions for SW and SHS fluids arc summarized in Section 



ITT1 The criterion followed in this paper to define the ef- 
fective stickiness parameter and some consequences are 
elaborated in Section ITTT1 Section llVl presents the com- 
parison between our MC simulations for SW and those 
reported in Ref. |^] for SHS, while the comparison with 
theoretical predictions is presented in Section The 
case of SW mixtures is considered in Section IvTl Finally, 
the paper ends with some concluding remarks. 



II. BASIC EQUATIONS 

In this Section we introduce the notation and display 
the basic equations that will be needed in the paper. 
A. The square-well fluid 

The SW interaction potential is given by Eq. (jX.ip . 
Henceforth, the distance r is assumed to be measured 
in units of the hard-core diameter ct, so that we take 
cr = l. The relevant physical information about the 
system is contained in the radial distribution function 
5sw( r l J 7j T*; A), where the notation emphasizes that it 
depends on the thermodynamic state (characterized by 
the packing fraction r\ and the reduced temperature T*) 
and on the relative range A. The corresponding cavity 
function is y S w(r\v, T*; A) = e vsv,( - r ^ kBT g sw {r \r), T*;X). 
Inverting this relationship, one has 



gsw(r\v,T*;X) 



0, r<l, 
e 1/T *2/sw(rkT*;A), 1< r < A, 
ysw(r\v,T*;X), r > A. 

Since the cavity function is continuous everywhere^ Eq. 
(|2.1() implies that the radial distribution is discontin- 
uous not only at r = 1, but also at r = A, namely 
ffsw (A- \n, T*; A)/ 5SW (A + |77, T*; A) — e 1 ^' . 

The structure factor is defined as S(k) = 1 + ph(k), 

where h(k) is the Fourier transform of the total correla- 
tion function h(r) = g(r) — 1. In the case of the SW fluid, 
S(k) can be obtained from y(r) as 



Ssw(k\i], T*;\) = ! + — J drrsmkr[g sw (r\r],T*;X)-l] 



= 1 -\ — — < — J dr rs'mkr + (e 1 ^ — 1) j dr rs'm kr 2/sw( r |??j T*; A) 



24ry 
~k 

"OO 

+ / drrsinfcr[y sw (r|77,T*;A)-l] \. (2.2) 
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The thermodynamic properties can also be expressed in terms of the cavity function. The virial equation yields 



Zsw(Vi T*\ A) = 1 



1 + 47/ 



4ry 



,1/T* 



dr r 3 g S w(r\v,T*; A) — <psw(r) 
1/Sw(% T*; A) — (e 1 ^* - l)X 3 y sw (X\ V: T*; A) 



(2.3) 



where Z = p/pksT is the compressibility factor, p being the pressure. The compressibility equation gives 
XswiVi T*;X) = S'sw(0|?7, T*; A) 

= 1 + 2477 1 (e 1/T * - 1) f drr 2 y sw (r\T 1 ,T*;\)+ I dr r 2 [ysw(r|r,, T*; A) - 1] - ~ 1 , (2.4) 



r 



where \ = ksT(dp/dp)T is the isothermal susceptibility. with the stickiness parameter ljl.21) kept constant: 
Finally, the excess internal energy per particle u is 



i-u sw (?7,T*; A) = ^ 



dr ipsw(r)gsw(r\r],T*; A) 



= -12r7e 1/T * / drr 2 y sw (r|77,T*;A). 

(2.5) 

Therefore, — 2uswA coincides with the coordination 
number. 



B. The sticky- hard-sphere limit 

The SHS system can be obtained from the SW one by 
formally taking the combined limits A — + 1 and T* — > 



USHs(r\r], t) = Jim ^ ysw(r\r), T*; A). (2.6) 



Equation (|2.1|l then reduces to 



9sns(r\r],T) = — 5(r - 1) + 6(r - 1) t/SHs(r|r?, r), 
izr 

so that the radial distribution function becomes identical 
to the cavity function for r > 1, but a delta-peak of am- 
plitude proportional to the stickiness appears at contact. 

In the SHS limit, Eqs. $F2&-$H& become 



Ssm(k\r],T) = 1 + ^ j- 



y drrsmfcr + y S Hs(l|»y,T) + y 



drrsin/cr [ysHs(f|?7, t) - 1] L (2.8) 



Zsm{v, T ) = 1 + 4t/ 



2/shs(%,t) - y^-ysHsUI 7 ?^) 



(2.9) 



f i z" 00 l 

XSHs(r?,r) = 1 + 24t? j— 2/shs(%,t) + J drr 2 [y S ns(r\v, r) - 1] - - 



(2.10) 



-Wshs(»7,t) = --2/SHs(l|?7, r )- (2-11) 



In Eq. H2.9fl . ?/(?") denotes the hrst derivative of y(r). 



III. EFFECTIVE "STICKINESS" 

Now we assume an SW fluid with a small, but non- 
zero, value of the well width A — 1. It is natural to ex- 
pect that its structural properties are, at least to some 
extent, close to those of an SHS fluid characterized by 
suitably chosen effective packing fraction r) e s and effec- 
tive stickiness t~q . Comparison between Eqs. (|2.1|) and 
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(|2.7() shows that the radial distribution function g(r) is 
not the adequate quantity to establish the approximate 
mapping SW^SHS. This also applies to the structure 
factor S(k) since it is essentially the Fourier counterpart 
of g(r). Therefore, we choose the cavity function y(r), 
which is a more regular function than g(r). Thus, we 
consider here the approximation 



the same as that of the SHS fluid, i.e., S| W (T*;A) 
Bf HS (T off (T*, A)). This yields 



Teff(T*,A) 



12(eVT' - l)a(A) 



(3.9) 



with a(A) = 03(A), where 



ysw( r \V: T*;X) w y S Hs{r\ri e s(v,T*;X),T cfi {r],T*;X)). 

(3.1) 

Note that this implies gsw(r\iliT*; A) 7^ 
9sns(r\Vcs{v,T*; A),T eff (?7, T*; A)) and Ssw(&|»7! r *5 A) 7^ 
-SsHS^keff^,? 1 *; A),T cff (?7,r*; A)). The next step 
consists of defining criteria to determine the effective 
quantities r] e s(ri, T*; A) and T e g(r), T*\ A)). There exist 
some obvious constraints. First, one must recover the 
hard-sphere (HS) cavity function from j/gw if cither 
A — > 1 at fixed T* or T* — > 00 at fixed A. More explicitly, 



lim y sw (r\r),T*;X) 



m.s(r\v), 



(3.2) 



a n (X) = -(A"-l). 
n 



(3.10) 



The form (|3.9|l with a(A) = a„(A) is consistent with the 
constraints l|3.5jl — l|3.7jl . regardless of the value of n. Since 
in the second virial coefficient one simply makes use of the 
zero-density cavity function y(r) — > 1, the choice a(X) = 
«3(A) is not necessarily the most efficient one from the 
point of view of Eq. (|3.1|) at finite density. Another com- 
mon choice^. 21 is r cff (T*, A) = e^T* /12(1 - A" 1 ), but 
it is inconsistent with the HS constraint (|3.t)|) . The same 
happens with the direct extrapolation of i|1.2|) at finite 
T*, i.e., r ofr (T*,A) = e- 1 / T 7l2(A-l). 

As a guide to choose r e ff(T*,A), let us resort to the 
exact cavity function to first order in density^ In that 
case one has 



lim ysw(r\rj,T*;X)=yns(r\r]), 



(3.3) 



ysw(r\ri, T*; A) = 1 + y^(r\T*;X)r, + 0( v 2 ), (3.11) 



lim ysm{r\r),T) = yns(r\v)- 

T — >00 

This implies that 



(3.4) 



lim r) eS (r), T*; A) = 77, lim T eff (rj, T*; A) = 00, (3.5) 



lim r] e s(r], T*; A) = 77, lim T e s(r], T*; X) = 00. 

T* — *oo T*— »oo 

(3.6) 

Moreover, Eq. Q3.1[l must become exact if the SHS limit 
is taken on the left-hand side and Eq. (|2.t)|) is applied. 
Thus, 

lim 12(A - l)e 1 / T *r cff (7 7 , T*; A) = 1. (3.7) 



ysns(r\v,r) = 1 + y^ s (r\r)T, + 0(r) 2 ), (3.12) 

where the expressions of ?/g^(r|T*;A) and 2/sHs( r l T ) 
are given in the Appendix. Comparison between both 
functions suggests to keep Eq. (|3.9|l . except that it 
remains to determine a(A) subject to the constraint 
limA— >i a(A)/(A — 1) = 1. In the region 1 < r < 2 we 
have 

y^{r\T*;X) = 8-6r+lr 3 

3(A 2 - l) 2 r- x - 8(A 3 - 1) + 6(A 2 - 1)7 



12r cff (T*,A)a(A) 
3(A 2 -l) 2 r-! 



[12r eff (T*,A)a(A)] 2 



(3.13) 



As is usualjS&iSiiSiijLSiiSS l e t us assume for simplicity 
that r] fi(r/,T*; A) = 7? e ff(??, A) is independent of temper- 
ature and r e ff (77, T*; X) — r e gf(T*;A) is independent of 
density. The first choice allows us to determine 77 e ff by 
applying the first equality of l|3.6|l , which holds regardless 
of the value of A > 1. It is then obvious that 



Vcsiv, A) = 7/. 



(3.8) 



On the other hand, the choice r] e s(r], A) = 77A 3 adopted 
by Menon et &\H is in conflict with the HS limit (|3.6() . 

In order to determine r e ff(T*,A), several criteria can 
be used. The conventional one 2ii2ii&2& is to impose 
that the second virial coefficient of the SW fluid be 



VshsWt) =8-6r+-r 3 



12r 2 



(3.14) 



It is clear that no unique choice of a(A) can make 

y { sw(r\T*;X) = y^ s (r|r off (T*, A)) for arbitrary r, T*, 
and A. In the limit of very low temperatures (T* — > 
=> T c ff — * 0) the respective last terms on the right- 
hand sides of Eqs. I|3.13|l and i|3.14|) dominate, so that 

Vsw( r \ T *' X ) 2/SHs( r l T eff( T *) A )) Provided th at on e 
takes a(A) = aa(A), where a„(A) is defined by Eq. I|3.10|l . 
Nevertheless, for temperatures such that 12r e g ~ 1, all 
the terms in Eqs. (|3.13|l and (|3.14|) are of the same order, 
so that again a(X) = 02(A) might not be the best choice. 
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FIG. 1: (Color online) Plot of the first-order virial coefficient 
of the cavity function, y^'(r), for an SW fluid with A = 1.05 
and T* = 0.35 and for SHS fluids with values of r given by 
Eqs. lOHJ and ll3~Ttjt with n = 3, n = 2, n = 1, n 0, and 
n = -1. 



Restricting ourselves, for the sake of simplicity, to func- 



tions a(A) of the form IJ3.10J) with n — integer, it turns 
out that the optimal choice is generally a(A) = ai(A) for 
temperatures close to the critical one, i.e., if r e g ~ 0.1,— 
This is illustrated in Fig. ^ for the case A = 1.05 and 
T* = 0.35. While the choice cv(A) = a 2 (A) = ±(A 2 - 1) 

gives a value of y^ s in excellent agreement with j/g^ 
at r = 1, the best global agreement is provided by the 
choice a(X) = ai(A) = A — 1. For larger temperatures 
the optimal choice could change to a(X) = ao(A) = In A 
or even a(X) = a_i(A) = 1 — A -1 . Since we want to 
keep a(X) independent of temperature, in what follows 
we adopt the choice a(A) = evi(A) for simplicity and also 
to enhance the cases with r c g rs 0.1. In summary, the 
definition of the effective stickiness parameter is taken as 

This is the same choice as considered by de Kruif et ali£. 
With this definition, the difference between j/ s ^ and j/g^s 
in the range 1 < r < 2 becomes 



y M(r\T*,\)-yg> s (r\T eS (T*,\)) = 



1 + 6(2 - 4r + r 2 )r eff 



(A-l) + 0((A-f) 2 ). 



(3.16) 



r 



Let us now explore the consequences of the approx- 
imation IJ3.1J1 . complemented by Eqs. (|3.8|l and [J3.15JI . 
on some of the physical properties. We first note that 
a simple approximate expression can be obtained for 
the excess internal energy of the SW fluid in terms of 
the value of the SHS cavity function and its slope at 
r = 1. First, we make the approximation ysw( r ) ~ 
2/Shs(1) + 2/sHs(l)( r — 1) inside the integral of Eq. IJ2.5JI . 
Next we perform the integral, make use of Eq. I|3.15|l and 
neglect terms nonlinear in A — 1. The result is 



1 



1 



Vmsi 1 ) 



-usw « — ushs 1+ 1 + 12r off + " SHS \* (A - 1) 
e e L V 2y S Hs(l)/ 

(3.17) 

where the arguments denoting the dependence of 
the quantities on 77, T* and A have been omitted 
for simplicity. Equation (|3.17() shows that, even if 
the approximation (J3.1J) is satisfied, usw(t)i T*; A) 7^ 
usm{v, T c s(T*;X)). 

Analogously, making yswO) sa 2/SHs(l)+2/s H s( 1 )( r_1 ) 
in the interval 1 < r < A, Eq. IJ2.2J) yields 



•S'sw(fc) - SsHs(fc) 



[ysHs(l)(A:cosfc + sin A:) 



kr e g 

+y SHS (l)sinfe](A-l), (3.18) 
where again terms nonlinear in A — 1 have been neglected 



and the arguments denoting the dependence of the quan- 
tities on 77, T* , and A have been omitted. Taking the limit 
k — > in the above expression, we have 



Xsw — XSHS 



V 

T e ff 



[2y SH s(l)+ysHs(l)](A-l). (3.19) 



Finally, from Eqs. (|2.3|l and (|2.9fl one gets 



Zsw — Zshs 



V 

TeS 



izshs(i) + TjysW 1 ) 



6' 



(A-l) 



(3.20) 

Equations H3.17JI - H3.20J) must be interpreted as simple 
heuristic approximations relating some physical proper- 
ties of SHS and short-range SW fluids. Thus, they do 
not give the first-order terms in a systematic expansion 
in powers of A — 1 because they are based on the ansatz 
(J3.1|) , which ignores terms of order A — 1 in the difference 
ysw( r ) — ysHs( r ), as shown by Eq. (J3.16I) . Note also that 
Eqs. (J3.17JI - IJ3.20JI depend on our choice IJ3.15JI for r e ff. 



IV. COMPARISON WITH SIMULATIONS OF 
STICKY HARD SPHERES 

We have performed conventional Monte Carlo simu- 
lations using the canonical ensemble and employing the 
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a The cases with A = 1 correspond to MC simulations for SHS 
performed by Miller and Frenkel 21 

TABLE I: SW systems simulated and corresponding values of 
2/(1), 2/'(l), and -u/e. 



linked cell-list method in a cubic box with a standard 
periodic boundary conditions. N = 4000 particles in- 
teracting via the SW potential l|l.lf) were displaced ac- 
cording to the Metropolis algorithm to create a sample 
of configurations. The calculation of the desired quanti- 
ties (internal energy and radial distribution function) was 
organized in cycles. Any randomly chosen particle was 
attempted to move and the cycle was repeated 100 x N 
times for equilibration before the system was analyzed. 
The process was repeated 1000 times to get one averaged 
set of values. Each run was divided into 20 such blocks, 
so that at least 10 9 configurations were generated. The 
distance of attempted displacements was a mixture of 
one value adjusted for each A so that the acceptance ra- 
tio was 10%-15% and another value ten times shorter to 
take care of tiny attractive effects. 

We have considered three different ranges: A = 1.05, 
1.02. and 1.01. In each case, a temperature T* has been 
chosen such that the effective stickiness parameter r Q g 
defined by Eq. (|oM5fl takes the same value. In order to 
test the approximation l|3.1(l , the values of r\ and r c g have 
been taken the same as those considered by Miller and 
Frenkel in Ref. 0. They are listed in Table [I] where 
also the internal energy —u/e and the values at contact 
y(l) and y'(l) are included. The two latter quantities 
have been obtained from a quadratic fit of the simulation 
data for ysw(0 and ysns{ r ) m the interval 1 < r < 1.1. 
Inasmuch as Eq. I|3.1|l is a reasonable approximation, one 
would expect that, at a given state point (77, T e fr), the four 
cavity functions, i.e., ysus{r) and ysw(f) with A = 1.05, 
1.02, and 1.01, overlap to some extent. Table [I] shows 
that the contact value of the cavity function, y(l), and 
its first derivative, ?/(l), are indeed quite similar for the 
four systems at each state. At the smallest temperature 
(r c ff = 0.13), the values of y(l) agree within statistical 
errors. At the intermediate temperature (t c s = 0.2), 
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FIG. 2: (Color online) Cavity function at 77 = 0.164 and r e a = 
0.13. The solid line corresponds to MC simulations for SHS 
by Miller and Frenkeli— The dotted, dashed, and dash-dotted 
curves correspond to our simulations for SW systems with 
A = 1.05, 1.02, and 1.01, respectively. 
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FIG. 3: (Color online) Same as in Fig.[T]but for 77 = 0.32 and 
t cS = 0.2. 
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FIG. 4: (Color online) Same as in Fig. 0but for 77 = 0.4 and 
T cft = 0.5. 
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2/sw(l) f° r ^ = 1-05 deviates from j/shs(1) less than 
1%, while ysw(l) f°r ^ = 1-02 and for A = 1.01 prac- 
tically coincide with ysHs(l)- At the highest tempera- 
ture (r c ff = 0.5), ysw(l) differs from j/shs(1) less than 
5%, 2%, and 1% for A = 1.05, 1.02, and 1.01, respec- 
tively. The fact that the agreement between ysw(l) and 
2/shs(1) slightly decreases as the temperature increases is 
associated with the choice of r e g given by Eq. (|3.15|) . As 
discussed in Sec. lIIII that choice was expected to be bet- 
ter for temperatures such that t c ^ ~ 0.1 than for larger 
ones. Similar conclusions can be drawn from the com- 
parison of 2/sw(l) ana - 2/shs(-1)) although in that case we 
have observed that the influence of the noise of the simu- 
lation data and that of the fitting procedure are slightly 
larger than in the case of y(l). 

Even if ysw{ r ) ~ ysns(f), that does not mean that 
it sw w ushs since the excess internal energy monitors 
the correlation function in the range 1 < r < A and so it 
is sensitive to the width of the well. In fact, Table[I]shows 
that the magnitude of the excess internal energy changes 
much more than the contact value of the cavity function 
when shrinking the well width. On the other hand, Eq. 
H3.17(l , which is implemented in Table [I] with the simu- 
lation data of J/shs(1) and ^shsW' incorporates terms 
linear in A — 1 and so it provides a very good estimate of 
the SW internal energy. 

A more complete comparison between ysw( r ) and 
2/SHs( r ) is provided by Figs. |2JB] There it is shown 
that the global shape of y sw (r|?7, T*; -V) at fixed 77 and 
t s{T* , A) is only weakly influenced by the value of 
A. In particular, the precursors of the singularities 
(delta-peaks and/or discontinuities )&2i of ysns( r ) at r = 

y/§, |, V3, 2, . . . are clearly apparent in the SW fluids, 
especially for A = 1.01. Apart from those local singu- 
larities, one can observe that the SW cavity function for 
A < 1.02 is reasonably well described by the SHS one, 
while slight but visible deviations are present in the case 
A = 1.05. 

It must be emphasized again that the good agreement 
between the SHS and SW cavity functions does not imply 
the same in the case of the radial distribution function 
inside the well (1 < r < A). As a consequence, the struc- 
ture factor, being related to the Fourier transform of h(r), 
can be expected to be affected by the difference between 
ffSHsM and gsw( r ) near r = 1. To visualize this effect, 
we compare in Figs. |3HZ|thc SHS and SW structure fac- 
tors, in the latter case with A = 1.02. Although both 
functions are very close each other in the three cases, 
some systematic differences are apparent. In particular, 
the maxima and minima of Ssw{k) are slightly more pro- 
nounced and shifted to the left than those of Ssns(k)- 
The difference AS(k) = S'sw(fc) — 5shs(&) is plotted in 
the insets of Figs.EHU Except for small k, it is observed 
that the simulation data of AS(k) are remarkably well 
described by the approximation 13.18fl . where the simu- 
lation values for ?/shs(1) and j/gnsW are used. This con- 
firms that the differences between ysw( r ) and ysus( r )i 




FIG. 5: (Color online) Structure factor S(k) at r\ = 0.164 and 
T c ff = 0.13. The solid line corresponds to MC simulations for 
SHS by Miller and Frenkeli-i The dashed curve corresponds 
to our simulations for an SW system with A = 1.02. The 
inset shows the difference AS(k) = Ssw(k) — Ssns(k) as ob- 
tained from simulations (solid line) and as given by Eq. (13.181 
(dashed line). 




FIG. 6: (Color online) Same as in Fig.|^]but for 77 = 0.32 and 
T cff = 0.2. 



which are neglected in Eq. I|3.18fl . have a practically ir- 
relevant influence on the difference between Ssw(k) and 

S'sHs(fc)- 



V. COMPARISON WITH THEORETICAL 
PREDICTIONS 

The approximate mapping ysw( r ) — * ysns( r ) can be 
exploited to get theoretical predictions for the radial dis- 
tribution function gsw( r ) of short-range SW fluids from 
the availability of analytical expressions for gsHsM- In 
particular, we will use Baxter's solution of the PY equa- 
tion for SHS^ as well as an extension of it proposed by 
two of usii^ For the sake of completeness, we will also 
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FIG. 7: (Color online) Same as in Fig. |S]but for r\ = 0.4 and 
T cft = 0.5. 



consider an alternative theoretical approach directly de- 
vised for SW fluids i^SiSL^S Before comparing with our 
MC simulations, the expressions stemming from these 
three theories are briefly described. 



A. Percus— Yevick and rational- function 
approximations for sticky hard spheres 



As is well known, Baxter was 
solution of the PY theory for the 
lution was later interpreted as the 
general class of approximations^ 
constructed in the Laplace space 
Laplace transform of rg(r): 

POO 

G(s) = / dr e~ rs rg(r) = s — 
Ji 1 

where the last equality defines 
F(s). In real space, 



able to find the exact 
case of SHS. 4 This so- 
simplest case of a more 
— Both approaches are 
. Let us introduce the 

F (s)e~ S f51) 
+ 12r ? F(s)e- s ' V ' ; 

the auxiliary function 



-. oo 



(5.2) 



n=l 



where f n (r) is the inverse Laplace transform of s[F(s)] n 
and O(x) is Heaviside's step function. The auxiliary func- 
tion F(s) must comply with some consistency conditions 
for small s and for large st^^ 



F(s) 



-1277 



0(s 5 ), 



F(s) = y(l)s-'(^ + s-^+0(s- 3 ). 



(5.3) 



(5.4) 



A simple class of approximations consists of assuming 
a rational form for F(s): 



F(s) 



(5.5) 



with m > 2. The simplest case corresponds to m = 2, 
the coefficients L%, L2, S%, S2, and S3 being explicitly 
determined 1 ^ from the application of the physical con- 
straints l|5. 3[) and (|5.4|) . This leads to the exact solution 
of the PY equation for SHSi^ One of the shortcomings 
of this PY-SHS theory is that it yields inconsistent equa- 
tions of state via the energy and compressibility routes. 

A more flexible approximation is obtained by setting 
m — 3 in Eq. (|5.5() and fixing the additional parameters 
L 3 and S4 by consistently imposing a given equation of 
state through the contact value y(l) and the isothermal 
susceptibility x~ We will refer to the ansatz l|5.5l) with 
L3 and 5*4 determined in this way as the rational-function 
approximation for SHS (RFA-SHS). Given a compress- 
ibility factor Zshs(?7,t) f° r the SHS fluid, the isother- 
mal compressibility and the excess internal energy are 
obtained as 1 ^ 

d 

XSHS = (VZSHS) , (5.6) 

usm d f v ,Z sm (ri',T) -1 

= -t— / dr] . (5.7) 

e dr J T}' 

The contact value ?/shs(1) consistent with Zshs can then 
be obtained from the internal energy through Eq. JSTHJ. 
For the compressibility factor .Zshs(?7, t) we will use the 
empirical form recently proposed by Miller and Frenkel. 22 



B. Rational- function approximations for 
square-well fluids 

In the case of SW fluids with a hard-core diameter a = 
1 and range A one can still define the auxiliary function 
F(s) trough Eq. H5.1|) . so that Eqs. (|5.2|l and l|5.3|l also 
hold in this casei 14 ^ However, the constraint (|5.4(l is 
now replaced by 

(5.8) 



F( S )=g(l + ) S -'+0(s-*). 



In the spirit of the ansatz (|5.5|l . the simplest rational- 
function approximation for SW systems (RFA-SW) has 
the form 2 S*2I 



F(s) 



1 1 + A + Kis- {A + K 2 s)e- {x -^ s 



12rj 



1 + Sis + S 2 s 2 + S 3 s 3 



(5.9) 



Application of the condition (|5.3|l allows one to express 
Ki, Si, S2, and ^3 as linear functions of A and Ki- Next, 
A is taken for simplicity as independent of density, so that 
A = e x l T — 1. Finally, K2 is obtained by enforcing the 
continuity of y(r) at r = A. 

It can be proved 2 ^ that in the SHS limit (A — > 1, T* — > 
with t = finite), one has ^(A— 1) — > L2 and ^1 — ^2 + 
A(A-l) — * Li, so that Eq. (|5.9|l reduces to Eq. 15.51) with 
m = 2, i.e., the PY-SHS solution. An improved RFA- 
SW theory that would reduce to the RFA-SHS theory 
would require the addition of new terms in the numerator 
and denominator of Eq. (|5.9|) , but then extra constraints 
would be needed to determine those new coefficients. 
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FIG. 8: (Color online) Radial distribution function of an SW 
fluid with A = 1.05 at r? = 0.164 and T* = 0.381. The 
circles correspond to our MC simulations and the solid lines 
represent the RFA-SW predictions. The results obtained from 
the approximation il-j. U with r e g = 0.13 by using the PY and 
RFA theories for the SHS fluid are represented by the dotted 
and dashed lines, respectively. The inset shows g(r) inside the 
well. Outside that region, the PY-SHS and RFA-SHS curves 
are practically indistinguishable. 
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FIG. 9: (Color online) Same as in Fig. H but for 77 = 0.32, 
T* = 0.448 and r cff = 0.2. 



FIG. 10: (Color online) Same as in Fig. H but for 77 = 0.4, 
T* = 0.682 and r cff = 0.5. 



empirical equation of state for SHS fluids^ as explained 
above. Despite the good theoretical behaviour, none of 
the theories predict the distortion around r = 1.77, which 
is the precursor of the singularity of the exact gsHS &t 
r = V3 » 1.73. Obviously, the PY-SHS and RFA-SHS 
approximations transfer to <7sw the discontinuity of gsns 
at r = 2, while the RFA-SW approximation accurately 
captures the rapid but continuous decay of <?shs near 
r = 2. In general, both SHS theories yield results which 
are practically indistinguishable, especially for the two 
highest temperatures (T* = 0.448 and T* = 0.682). Only 
inside the attractive well (1 < r < 1.05) can one distin- 
guish the PY-SHS and RFA-SHS predictions, as shown in 
the insets. This region is especially important for SW flu- 
ids since it is directly related to the coordination number 
and the excess internal energy. We observe that the best 
performance in the region 1 < r < 1.05 is done by the 
RFA-SHS theory at the lowest temperature (T* = 0.381) 
and by the RFA-SW theory at the other two tempera- 
tures (T* = 0.448 and T* = 0.682). In thi s respect, it 
must be borne in mind that the choice (|3.15|) of the effec- 
tive stickiness parameter is more appropriate for temper- 
atures such that T e s ~ 0.1 than for higher temperatures. 



C. Comparison 

We are now in conditions of comparing the three the- 
ories against our MC simulation data. We restrict our- 
selves to the case A = 1.05. Since the radial distribu- 
tion function g(r) is a more direct and intuitive quantity 
than the cavity function y(r), we present the plots in 
terms of g(r) in Figs. ISHTUI We observe a good general 
agreement between theories and simulations for the three 
states. This is noteworthy since no fitting has been em- 
ployed. Only in the case of the RFA-SHS approach there 
exist two free parameters (L3 and S4) which are deter- 
mined by enforcing thermodynamic consistency with an 



VI. MIXTURES OF SW FLUIDS 

Let us now consider the case of multicomponent sys- 
tems of particles interacting via SW potentials: 

!oo, r < crij, 
-Eij, Oij < r < XijCTij, (6.1) 
0, r>Xij(Tij. 

In the one-component case, there are three parameters 
characterizing the interaction (er, e, and A) and two pa- 
rameters defining the thermodynamic state (p and T). 
Without loss of generality, we can use a and e to fix the 
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units of distance and energy, respectively, so that only 
three independent dimensionless parameters remain (A, 
77, and T*). In the case of mixtures, however, the number 
of parameters is in general much higher. Restricting our- 
selves to the binary case, we have three thermodynamic 
quantities (the total density p, the mole fraction x\ of one 
of the species, and the temperature T) plus nine interac- 
tion parameters (the diameters cr^ , the well depths , 
and the relative ranges Ay). One of the diameters can 
be used to define the length unit and one of the depths 
can be used to define the energy unit, so that the number 
of independent parameters defining the problem is ten in 
general^ 

The SW potentials (|f> . lft become SHS potentials if 
Ay — > 1 and eij — > 00 by keeping constant the parameters 



« 12(e e -/^ T -l)(Ay-l)' K ' 
I 



Class 


(Til 


0"22 


0"12 


An 


A22 


A12 


en 


£22 


ei2 xi 


I 

II 


1 
1 


1 

3 


1 

2 


1 

A 


1 

A 


A 
A 




e 




c 


- \ 



TABLE II: Classes of SW mixtures considered in the simula- 
tions. 



where yij(r) = e v ^ ( - r ^ kBT g ij (r). 

Given the high dimensionality of the parameter space 
for binary SW mixtures, we restrict ourselves to the two 
classes of mixtures shown in Table ITT1 Class I corresponds 
to equimolar and symmetric mixtures in which the self 
interactions 1-1 and 2-2 are of HS type, while the cross 
interaction 1-2 is of SW type. Therefore, rff = r|f = 00 
and t°2 — T cS- This class of mixtures can be used to de- 
scribe adsorption phenomenal Class II is representative 
of asymmetric hard-core diameters, but with common 
well depths and relative well widths, so that r,? ff = T e s- 
It only remains to fix for each class the values of the 
range A, the reduced temperature T* — ksT/e and the 
reduced density pa 3 or, equivalently, the total packing 
fraction i] = np(xi<jfi + x-ioW) /6. The values considered 
in this paper for those quantities are shown in Table ITTT1 
This yields a total of 12 different mixtures simulated. For 
each A, the temperature is chosen so that r e g = 0.2 and 
0.5, in analogy with some of the cases considered in Table 

□ 

Now we compare, for each class and each value of 
T e ff, the cavity functions yij(r) obtained from our MC 
simulations. Given the scarcity of simulations for SHS 
mixtures we use the analytical solution of the PY 



Now, each pair (ey,Ay) gives birth to a parameter r° ff 
and so the general number of parameters characterizing 
the binary mixture reduces from ten to seven. Note that, 
although the term —1 appearing in the denominator of 
Eq. (|6.2|l can be neglected in the SHS limit, we have kept 
it so that Eq. (|6.2|l defines the generalization to mixtures 
of the effective stickiness parameter defined by Eq. (f 3 . 1 f>|> . 
The exact solution of the PY equation for SHS mixtures 
is knownii in the additive case cry = (an + (Xjj)/2. Sys- 
tems of SHS mixtures have been considered by several 
authorsi 7 ! 11 ! 17 ! 18 ' 23 ! 24 An extension of the RFA described 
in Subscc. rV Al to SHS mixtures is also available^ 



In agreement with the philosophy behind the approxi- 
mation (|3.1|) for the one-component case, we can expect 
that, for narrow SW potentials, 



(6.3) 

I 



_n Tcff a t* 

0.4 0.2 1.05 0.448 

1.02 0.324 

1.01 0.266 
0.5 1.05 0.682 

1.02 0.448 
1.01 0.348 



TABLE III: Values of 77, A and T* considered in the simula- 
tions. 



closured We cannot use the expressions provided by the 
RFA method 3 ^ because the latter relies upon the knowl- 
edge of the equation of state for SHS mixtures. While 
an accurate empirical equation of state exists for one- 
component SHS fluids^ its extension to mixtures is not 
known. 

The two states corresponding to mixtures of class I are 
considered in Figs.llllandll2l We observe that the range 
A = 1.05 is again not small enough to produce a good 
overlap with the cases A = 1.02 and A = 1.01. These 
two latter SW mixtures can be expected to be close to 
the true SHS mixture. In this respect, the deviations of 
our SW simulations with A = 1.02 and A = 1.01 from the 
SHS-PY theory reveal the limitations of the latter i-4 Fig- 
ures El an d El show that the PY theory describes better 
the pair correlation function yu = 7/22 of the particles 
interacting via the HS potential than the pair correla- 
tion function yi2 of the particles feeling a short-range 
attraction. Accordingly, the agreement improves when 



Vij {r\p,xi,T,{a k ^},{eu},{Xu}) s» y 4j (r\p, x x , {a ki }, {r^ }), 
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r 

FIG. 11: (Color online) Cavity functions at r\ = 0.4 and r e a = 
0.2 for a binary mixture of class I (see Table |nj. The solid 
lines correspond to the analytical solution of the PY theory for 
SHS mixtures. The dotted, dashed, and dash-dotted curves 
correspond to our simulations for SW systems with A = 1.05, 
1.02, and 1.01, respectively. 



the temperature increases (from r c g = 0.2 to r c ff = 0.5) 
and hence the stickiness decreases. 

The results corresponding to the asymmetric mixtures 
of class II are plotted in Figs.ll3landll4l The convergence 
SW^SHS is similar to that in the symmetric cases since 
the deviations between the curves corresponding to the 
ranges A = 1.02 and A = 1.01 are again small. However, 
the performance of the PY theory is worse in Figs. El 
and El than in Figs. ^2 and El This could have been 
expected by the following arguments. On the one hand, 
the presence of stickiness is larger in class II than in class 
I. On the other hand, it is known that the PY theory is 
more accurate for one-component HS fluids (class I in the 
limit r — > oo) than for HS mixtures (class II in the limit 



VII. SUMMARY AND CONCLUDING 
REMARKS 

In this paper we have investigated the possibility of 
representing the structural and thermodynamic proper- 
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r 

FIG. 12: (Color online) Same as in Fig. El but for r? = 0.4 
and T e ff = 0.5. 



ties of short-range SW fluids by those of SHS fluids. It 
is sometimes assumed in the literature on colloidal sus- 
pensions that a range A < 1.1 is small enough to replace 
the more realistic SW model by the simpler SHS model. 
Moreover, the mapping SW— *SHS is usually assumed to 
hold at the level of the structure factor, i.e., Ssw(k) ~ 
Ssns(k). This in turn implies g S w{r) ~ SshsM, what 
misses the fact that gsw(r) has a jump discontinuity at 
r = A which becomes a delta singularity of <?SHs( r ) a t 
r = l+. 

Here, however, the ansatz Ssw(k) « S'sHs(fc) or , 
equivalently, gsw( r ) ~ 3shs(?") has been replaced by 
Usw{ r ) ~ ysns( r )- This still leaves open the problem 
of finding the effective packing fraction r) e g and the effec- 
tive stickiness parameter r c g of the SHS fluid that best 
mimics the physical properties of an SW fluid of range A 
at a given state (rj,T*). For simplicity, we have chosen 
r) e g to be independent of T* and T e g to be independent of 
T}. Next, the constraint of recovering the HS fluid in the 
limit T* — > oo leads to r/ c g — rj in a natural way. Since 
r c ff is assumed to be independent of 77, we have taken 
advantage of the exact knowledge of ysw{r) and j/SHs( r ) 
to first order in density as a guide to choose T e s- Con- 
sidering several possibilities of the form l|3.9|l . we have 
found that the optimal choice for temperatures such that 
r e ff ~ 0.1 is provided by Eq. Q3.15)) . as illustrated by Fig. 
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FIG. 13: (Color online) Cavity functions at r\ = 0.4 and T e g — 
0.2 for a binary mixture of class II (see Table HTl . The solid 
lines correspond to the analytical solution of the PY theory for 
SHS mixtures. The dotted, dashed, and dash-dotted curves 
correspond to our simulations for SW systems with A = 1.05, 
1.02, and 1.01, respectively. 



^ This choice differs from the conventional nej2*i2*i2i2£ 

based on the equality of the SW and SHS second virial 
coefficients, i.e., Bf w (T*;A) = Sf HS (r cf f (T*, A)). 

The condition i?| w = Bf HS is directly related to the 
assumption that the thermodynamic properties of the 
SW and SHS fluids are the same, at least for low den- 
sities. Our approach differs from this one. We have as- 
sumed that the cavity function is approximately the same 
in both systems and from this ansatz we have obtained 



FIG. 14: (Color online) Same as in Fig. EU but for r? = 0.4 
and T c ff = 0.5. 



expressions for the structure factor, the internal energy, 
the isothermal compressibility, and the pressure of the 
SW fluid in terms of quantities related to the SHS fluid. 
Those approximate relations have been checked by sim- 
ulation data in the case of the structure factor (see Figs. 
|SH7|) and the internal energy (see Table P). 

We have performed MC simulations of one-component 
SW fluids for A = 1.05, 1.02, and 1.01. For each case, 
we have considered the densities and temperatures indi- 
cated in Table [I] They correspond to the same values 
of the packing fraction and effective stickiness as consid- 
ered by Miller and Frenkel in their simulations of SHS 
fluidsiSi The comparisons (see Figs. show that the 
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SW cavity functions with A = 1.05 present small but vis- 
ible differences with respect to the SHS ones. On the 
other hand, the SW curves corresponding to A = 1.02 
and 1.01 have almost collapsed to the SHS curves, ex- 
hibiting the precursors of the singularities of ysns( r ) a t 

r = \J~^i §j V^3j 2, . . . 

While the PY equation has an exact solution for SHS, 
its solution for SW is not known exactly. However, we 
have exploited the knowledge of the PY-SHS solution, 4 
as well as of the rational-function approximation (RFA) 
for SHS, 15 to estimate the radial distribution function 
jsw(r) based on the ansatz j/sHs( r ) ~ ysw(r)- The re- 
sults for A — 1.05 show a good general agreement with 
the MC data (see Figs. l8* HTUj) . except near r — 2, where 
gsw(r) decays rapidly but does not present the jump dis- 
continuity of gsns( r )- From that point of view, a better 
agreement is provided by the RFA-SW theory^ 

In order to offer a wider perspective, we have also 
considered two classes of binary mixtures (see Table [HJl . 
Class I defines a symmetric equimolar mixture where the 
interactions among particles of the same species are of 
HS type, while those among particles of different species 
are of SW type. Mixtures of class II are asymmetric 
with a hard-core ratio 1:3 but with an attractive interac- 
tion of common depth and (relative) width. Comparison 
between our simulation data for SW mixtures and the 
theoretical PY-SHS predictions^ shows that the latter 
behaves better for mixtures of class I (see Figs. ^2 and 
I12|l than for mixtures of class II (see Figs.ll3land ll4fl . In 
each case, the agreement improves as the temperature in- 
creases (and hence the stickiness decreases). Apart from 
this comparison with the PY-SHS theory, our MC simu- 
lations show again that the range A = 1.05 is not small 
enough to get a good collapse with the cases A = 1.02 
and 1.01. 

One of the open avenues not fully explored in this pa- 
per refers to the analysis of the equation of state for short- 



range SW fluids obtained from the empirical equation 
of state for SHS fluids recently proposed by Miller and 
FrenkeliSa Equations (|3~T7jl . l|3~T3|l . and (j?T2T7)) provide 
three (approximate) alternative routes to the SW equa- 
tion of state. They require the knowledge of 2/shs(1)j 
2/g HS (l), and j/gHsW- The two former quantities can be 
obtained from the SHS equation of state through Eqs. 
(|2.9|) and (|2.11|) . while 2/shs(1) can be obtained from the 
RFA-SHS theory^ We plan to undertake this study in 
the near future. 
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APPENDIX: THE CAVITY FUNCTION TO 
FIRST ORDER IN DENSITY 

To first order in density, the cavity function of the SW 
fluid is given by Eq. (|3.11|l . The coefficient ia~ 



Vsw(r\T*;\) = A(r) + (e 1 ^* - l)B S w(r|A) + (e 1/T * - l) 2 C SW (r|A), (A.l) 

where (for A < 3) 

m ^ [8-6r+ ±r 3 , r < 2, 

Air) = < 2 (A. 2) 

10, r>2, 



BswHA) 



-12r + r 3 , 
3(A 2 - l) 2 r-! - 8(A 3 - 1) + 6(A 2 - l)r, 



3(A 2 
0, 



Ifr- 1 - 8(A 3 + 1) + 6(A 2 + l)r - r 3 



r < A- 1, 
A - 1 < r < 2, 
2 < r < A + 1, 
r > A + l, 



(A.3) 



14 



C sw (r\X) = < 



8(A 3 - 1) -6(A 2 + l)r + 
3(A 2 ~l) 2 r~\ 
3(A 2 -l) 2 r- 1 -8 
8A 3 -6A 2 r + ±r 3 , 



6r - ±r 3 



r < A — 1, 
A - 1 < r < 2, 
2 < r < A + 1, 
A + 1< r < 2A, 
r > 2A. 



(A.4) 



Taking the SHS limit, Eq. IjA.ljl reduces to 



yHs(r\r)=A(r) + — Bsus(r) + 7T ^C S u S (r), (A.5) 



where 



(12r) 2 



Jl2r-\ 0<r<2, 
CsHs(r) = 0, r > 2. 



(A.7) 



s SHS «= ; 2(r - 2) ' °;r 2 ' ^ 

0, r > 2, 
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